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The Woods-Saxon basis has been suggested to replace the widely used harmonic oscillator basis for 
solving the relativistic mean field (RMF) theory in order to generalize it to study exotic nuclei. As 
examples, relativistic Hartree theory is solved for spherical nuclei in a Woods-Saxon basis obtained by 
solving either the Schrodinger equation or the Dirac equation (labelled as SRHSWS and SRHDWS, 
respectively and SRHWS for both). In SRHDWS, the negative levels in the Dirac Sea must be 
properly included. The basis in SRHDWS could be smaller than that in SRHSWS which will 
simplify the deformed problem. The results from SRHWS are compared in detail with those from 
solving the spherical relativistic Hartree theory in the harmonic oscillator basis (SRHHO) and those 
in the coordinate space (SRHR). All of these approaches give identical nuclear properties such as 
total binding energies and root mean square radii for stable nuclei. For exotic nuclei, e.g., 72 Ca, 
SRHWS satisfactorily reproduces the neutron density distribution from SRHR, while SRHHO fails. 
It is shown that the Woods-Saxon basis can be extended to more complicated situations for exotic 
nuclei where both deformation and pairing have to be taken into account. 
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I. INTRODUCTION 

The existence of an average field in atomic nuclei re- 
vealed by the exceptional role of the nuclear magic num- 
bers provides the foundation of the nuclear shell model 
and various mean field approaches Q, 0, • This aver- 
age field is believed to be approximated most closely by 
a Woods-Saxon (WS) potential [4( either from analyzing 
the radial dependence of the nuclear force or by deriving 
it from a microscopic two-body force. 

Since the eigenf unctions for the WS potential can not 
be given analytically, as good approximations for stable 
nuclei, one often adopts the harmonic oscillator (HO) 
potential or the square well, in particular the former, 
in shell model calculations for both spherical |l| and de- 
formed nuclei The HO eigenfunctions also often serve 
as a complete basis in solving equations in both non- 
relativistic and relativistic mean field approximations, 
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such as the Skyrme Hartree- Fock (SHF), Hartree- Fock- 
Bogoliubov (HFB), relativistic Hartree (RH) and rela- 
tivistic Hartree-Bogoliubov (RHB) theories. In these ap- 
proaches, solution of corresponding equations is trans- 
formed to a matrix diagnolization problem which can be 
easily dealt with. 

However, due to the incorrect asymptotic property of 
the HO wave functions, the expansion in the localized 
HO basis are not appropriate for the description of drip 
line nuclei 0, 0, @ which display many interesting fea- 
tures because of the extremely weakly bound property, 
e.g., the coupling between bound states and the contin- 
uum due to the pairing correlation, large spacial density 
distributions, possible modifications of shell structure, el 
al. One must improve the asymptotic behavior of HO 
wave functions, e.g., by performing a local scaling trans- 
formation |9j- However, one does not know the scaling 
parameter beforehand thus the predictive power in this 
method is lost. 

A proper representation to solve the HFB or RHB 
equations for drip line nuclei is the coordinate space 0, 
[XOL [ill ll^ where wave functions are approximated on a 
spatial lattice and the continuum is discretized by suit- 
ably large box boundary conditions. The HFB method 
solved in r space can take fully into account all the mean- 
field effects of the coupling to the continuum [E 01 El El- 
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Nevertheless for deformed nuclei, working in r space be- 
comes much more difficult and numerically very sophis- 
ticated @- Particularly, it become very time consuming 
when the pairing correlation is included. Therefore much 
effort is made towards a more efficient solution of HFB or 
RHB equations, e.g., using natural orbitals 0] or work- 
ing on basis-spline Galerkin lattices 0], etc.. 

A reconciler between the HO basis and the r space 
may be the WS basis because (i) the WS potential repre- 
sents the nuclear average field more suitably than the HO 
potential and (ii) in principal there is no localization re- 
strictions in the WS potential. Although analytical wave 
functions can not be given for the WS potential, one may 
easily find numerical solutions for a spherical WS poten- 
tial in the r space by virtue of various effective methods 
of solving ordinary differential equations |l6l |. One can 
still use a large box boundary condition to discretize the 
continuum. These WS wave functions can thus be used 
as a complete basis for spherical or deformed systems and 
one finally comes back to the familiar matrix diagnoliza- 
tion problem. 

In the present work we restrict the application of this 
method to nuclei with spherical symmetry which largely 
facilitates the discussion of basic principles and allows 
presenting illustrations for the radial dependence of all 
relevant physical quantities like density distributions. We 
combine this approach with the relativistic Hartree the- 
ory 01 which provides a framework for describing the 
nuclear many body problem as a relativistic system of 
baryons and mesons and, together with its extensions 
with deformation and/or pairing included, have been suc- 
cessfully applied in calculations of nuclear matter and 
pro pert ies of finite nuclei throughout the periodic ta- 
ble llEO- 

The paper is organized as follows. In Sec. [HJ we give 
a brief reminder of the formalism of relativistic Hartree 
theory. The numerical details of solving it in the WS 
basis are given in Sec. II I II In Sec. IIVI we present our 
results and compare them with those obtained in the HO 
basis and in the r space. We also discuss the contribution 
from negative levels in the Dirac sea in the same section. 
Finally, the work is summarized in Sec. 

Throughout the paper, the relativistic Hartree theo- 
ries solved in the r space, in the HO basis and in the 
WS basis are abbreviated as "SRHR" , "SRHHO" and 
"SRHWS" where the first "S" represents "spherical". 
We use "SWS" and "DWS" to distinguish the WS ba- 
sis which is obtained from solving the Schrodinger equa- 
tion or the Dirac equation with initial WS potentials, 
respectively. Thus we have "SRHSWS" and "SRHDWS" 
theories. 



II. BASIC FORMALISM OF RELATIVISTIC 
HARTREE THEORY 

The starting point of the relativistic Hartree theory 
is a Lagrangian density where nucleons are described as 



Dirac spinors which interact via the exchange of several 
mesons (a, u, and p) and the photon [l7Lll8l lT9|. 
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with the summation convention used and the summation 
over i runs over all nucleons, j: = 7 M x M = l^x^, M the 
nucleon mass, and m a , g ai m w , <? w , m p , g p masses and 
coupling constants of the respective mesons. The non- 
linear self-coupling for the scalar mesons is given by [20l | 
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U(<t) = -m a a + ycr J + —a , 



(2) 



and field tensors for the vector mesons and the photon 
fields are defined as 



QfM/ — d^LO u — d v U[j,, 

R^v = d^pv - d v p^ - g p {p p x p v ), 
F pu — d^A v — d v A p . 



(3) 



The classical variation principle gives the equations of 
motion for the nucleons, mesons and the photon. As in 
many applications, we study the ground state properties 
of nuclei with time reversal symmetry, thus the nucleon 
spinors are the eigenvectors of the stationary Dirac equa- 
tion 

[a ■ p + V(r) + (3(M + S(r))} ^(r) = e^i(r), (4) 

and equations of motion for the mesons and the photon 
are 

{-A + d a U{a))a{r) = -g*p a (r), 
(-A + mi)u°(r) = g^ Pv (r), 
(-A + m 2 p ) p°{r) = g P p 3 (r), 
-AA°(r) = ep p (r), 



(5) 



where uj q and A are time-like components of the vec- 
tor uj and the photon fields and p° the 3-component of 
the time-like component of the iso-vector vector p me- 
son. Equations and JSJ are coupled by the vector and 
scalar potentials 



J V(r) = gul oj°{r) +g p T 3 p°(r) + e 
\ S(r) = g a a(r), 

and various densities 
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For spherical nuclei, meson fields and densities depend 
only on the radial coordinate r. The spinor is character- 
ized by the angular momentum quantum numbers 
to, the parity, the isospin t = ±1/2 ("+" for neutrons 
and "— " for protons), and the radial quantum number 
a. The Dirac spinor has the form 



F2(r) 



Y, 



Xt a (t), j = l±~, 

(8) 

with G*(r)/r and F£(r)/r the radial wave functions for 
the upper and lower components and Yj m (9,4>) the spin 
spherical harmonics where k — (— l)J+ i + 1 / 2 (j + 1/2) and 
1 = 1 + (— 1) J+ ' -1 ' 2 . The value of k of the upper compo- 
nent is used to label a state both for normal levels in the 
Fermi sea and for negative ones in the Dirac sea. States 
with the same k form a "block". The radial equation of 
the Dirac spinor, Eq. is reduced as 



TABLE I: Dependence of the average single particle energy, 
root mean square (rms) radius, and (r 4 ) 1 ' 4 on the mesh size 
Ar for the SRHSWS theory. The meson and Coulomb fields 
are obtained from SRHR calculations with the parameter set 
NL3, Ar = 0.05 fm, 7? max = 30 fm for 16 and 35 fm for 
208 Pb. In SRHSWS calculations, the parameter set NL3 is 
used. For 16 0, R msx = 4r A 1/3 = 12.8 fm, E cut = 300 MeV. 
For 208 Pb, ,R max = 3r A 1/3 = 22.6 fm, E cut = 200 MeV. The 
first row gives results from solving the Dirac equation in the 
coordinate space. 
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The meson field equations become simply radial Laplace 
equations of the form 

TTitfi are the meson masses for <fi = a, u>, p and zero for the 
photon. The source terms are 



s^O) 
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with 
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The above coupled equations have been solved in the r 
space [2l| and in the HO basis with the no sea and the 
mean field approximation. Here we depict briefly the pro- 
cedure of solving these coupled equations. With a set of 
estimated meson and photon fields, the scalar and vector 
potentials are calculated and the radial Dirac equation is 
solved. Thus obtained nucleon wave functions are used 
to calculate the source term of each radial Laplace equa- 
tion for mesons and the photon. New meson and photon 
fields are calculated from solving these Laplace equations. 
This procedure is iterated until a demanded accuracy 
is achieved. Laplace equations are usually solved using 



the Green's function method |21j,|22j though in Ref. |22| 
Laplace equations for mesons are solved in the HO ba- 
sis. SRHR, SRHHO and SRHWS differ from each other 
mainly in how to solve the Dirac equation. In the follow- 
ing the numerical solution of the Dirac equation in the 
WS basis will be presented. 



III. SOLVING THE DIRAC EQUATION IN A 
WOODS-SAXON BASIS AND NUMERICAL 
DETAILS 

A. Woods-Saxon basis from solving a Schrodinger 
equation (the SWS basis) 

For the Schrodinger equation with a spherical Woods- 
Saxon potential 



To 



Vws(r) 



I _|_ g(r— Ro)/ao ' 
OO, 



(13) 



where i? max is introduced for practical reasons to define 
the box boundary, the eigenfunction can be written as 
4>nlmi(r) = Rni(r)Yi mi (9,(f>). Its radial Schrodinger equa- 
tion is derived as 



1 9 2 d 
■r 



~2M \r 2 i),-' Or 
E n iR nl (r). 



1(1 + 1) 



Vws (r) 



Rni{r) 
(14) 



Equation (|14|) is solved on a discretized radial mesh 
with a mesh size Ar. R max (Ar) should be chosen larger 
(smaller) enough to make sure that the final results do 
not depend on it. The radial wave functions thus ob- 
tained form a complete basis, 



{R„i{r);n = 0, 1, • • • ; I = 0, 1, • • • ,n} 



(15) 



TABLE II: Dependence of the average single particle energy, 
rms radius, and {r 4 ) 1//4 on the difference An = n max — n max 
for the SRHSWS theory. The meson and Coulomb fields 
are obtained from SRHR calculations with the parameter set 
NL3, Ar = 0.05 fm, R max = 30 fm for ls O and 35 fm for 
208 Pb. In SRHSWS calculations, the parameter set NL3 is 
used. For ie O, i? max = 4r A 1/3 = 12.8 fm. For 208 Pb, 
i?max = 3r A 1/3 = 22.6 fm. 
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in terms of which the radial part of the upper and the 
lower components of the Dirac spinor in Eq. @ are ex- 
panded respectively as 
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The radial Dirac equation, Eq. JUJ, is transformed into 
the WS basis as 
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where the matrix elements are calculated as follows 

flmax 
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r 2 drR^(r) (V(r) - S(r) - M) R hl {r). 

In practical calculations, an energy cutoff E cut (relative 
to the nucleon mass M) is used to determine the cutoff 
of the radial quantum number n max for each block. In 
the expansion of the corresponding lower component, we 
take n max 

spurious states |22j. 

The following Woods-Saxon parameters have been 
used according to Ref. 



n™! An with An > 1 in order to avoid 



Vq = (-51 ± 33(iV - Z)/A) MeV, 
Ro = 1.27A 1 / 3 fm, a = 0.67 fm, 



(19) 
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FIG. 1: Deviations of the average single particle energy 
E ap /A (the upper panel), rms radius (r) 1//2 (the middle 
panel), and (r 4 ) 1 / 4 (the lower panel) of 16 O from the standard 
results versus the cutoff energy _E cut with different box size 
Rmax for the SRHSWS theory. The meson and Coulomb fields 
are obtained from SRHR calculations with the parameter set 
NL3, Ar = 0.05 fm, i? max = 30 fm. In SRHSWS calculations, 
the parameter set NL3 is used. For E cut = 100, 200, 300, 400 
MeV, the approximate maximum principal quantum number 
in each basis, iV max = 2n max + /, is given on the top of the 
plot. 



where '+' is for the neutron and ' — ' for the proton. As ex- 
pected, the dependence of final results on the initial WS 
potential is almost negligible. For example, a variation 
of Vb by 50% gives differences in total binding energies 
by less than 0.1% and charge radii by less than 0.5% for 
16 0, 48 Ca and 208 Pb. Such situation is also checked to be 
true for the other two parameters in the WS potential, 
Rq and oq. 



B. Woods-Saxon basis from solving a Dirac 
equation (the DWS basis) 

The radial Dirac equation, Eq. may be solved in 
the r space |21j with Woods-Saxon-like potentials for 
Vo (r) ± So (r) |23| within a spherical box of the size i? ma x, 
together with the spherical spinor which gives a complete 
WS basis 
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FIG. 2: Deviations of the average single particle energy 
E sp /A (the upper panel), rms radius (r) 1//2 (the middle 
panel), and (r 4 ) 1 ^ 4 (the lower panel) of 208 Pb from the stan- 
dard results versus the cutoff energy _E cut with different box 
size i?max for the SRHSWS theory. The meson and Coulomb 
fields are obtained from SRHR calculations with the parame- 
ter set NL3, Ar = 0.05 fm, i? max = 35 fm. In SRHSWS calcu- 
lations, the parameter set NL3 is used. For E cut = 100, 200, 
300, 400 MeV, the approximate maximum principal quantum 
number in each basis, iV max = 2n max + 1, is given on the top 
of the plot. 



FIG. 3: Deviations of the average single particle energy 
E ap /A (the upper panel), rms radius (r) 1//2 (the middle 
panel), and (r 4 ) 1 / 4 (the lower panel) of 16 from the stan- 
dard results versus the cutoff energy _E cut with different box 
size -R m ax for the SRHDWS theory. The meson and Coulomb 
fields are obtained from SRHR calculations with the parame- 
ter set NL3, Ar = 0.05 fm, i? max = 30 fm. In SRHDWS cal- 
culations, the parameter set NL3 is used. For E cut = 100, 200, 
300, 400 MeV, the approximate maximum principal quantum 
number in each basis, iVmax = 2n max + I, is given on the top 
of the plot. 



with n = 0, 1, • • • , K = ±1, ±2, • ■ • , and m — —j K , • • ■ , j K . 
■0° Km (r, s, i) takes the form of Eq. ©. We note that 
states both in the Fermi sea and in the Dirac sea should 
be included in the basis for the completeness. The nu- 
cleon wave function, Eq. (|SJ), can be expanded in terms 
of this set of basis as 

™max 

ip aKm (r,s,t) = c^Lm^M), (21) 

n=0 

where n max = nm ax +« max + l and the summation is over 
normal levels in the Fermi sea for < n < n^ ax and over 
negative levels in the Dirac sea for n+ ax + 1 < n < 
The negative states is obtained with the same method 
as the positive ones [2l]. In this WS basis, the Dirac 
equation, Eq. turns out to be 

^max 

Cajn^rri "T" ^ ^an^rnn — ^■cx^otmi ^ — E ' ' ' ) ^maxi (22) 
n=0 



with 

H' mn = «(r)| [AV(r) + (3AS(r)} |^(r)) 

= / " drGl (r) [Ay (r) + AS(r)} G° n (r) 
Jo 

+ f " drF^(r)[AV(r) - AS(r)] F°(r),(23) 
Jo 

where AV(r) = V(r)-V (r) and AS(r) = S(r)-S (r). 
The angular, spin, and isospin quantum numbers are 
omitted for brevity. 

It should be mentioned that Eq. © can be solved di- 
rectly in the r space with the same method of generating 
the DWS basis. It is our aim to test the validity of an 
efficient solution not only for the spherical RH model but 
also for its extension to include the deformation and/or 
the pairing correlation. In fact, if only the SRH theory 
is concerned, this procedure is just a replacement of the 
direct solution in the r space by a diagnolization of a 
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FIG. 4: Deviations of the average single particle energy 
E sp /A (the upper panel), rms radius (r) 1//2 (the middle 



panel), and (: 



Ml/4 



(the lower panel) of Pb from the stan- 



dard results versus the cutoff energy _E cut with different box 
size -R m ax for the SRHDWS theory. The meson and Coulomb 
fields are obtained from SRHR calculations with the parame- 
ter set NL3, Ar = 0.05 fm, i? max = 35 fm. In SRHDWS cal- 
culations, the parameter set NL3 is used. For E cut = 100, 200, 
300, 400 MeV, the approximate maximum principal quantum 
number in each basis, iV max = 2n max + I, is given on the top 
of the plot. 



matrix with some complication introduced by the fact 
that contributions from states in the Dirac sea must be 
included. 

An energy cutoff -E C ut (relative to the nucleon mass 
M) and the cutoff of radial quantum numbers n max are 
applied to normal levels alternatively according to prac- 
tical convenience. For the initial Woods-Saxon potentials 
V (r) ± S (r), we follow Ref. 0. 



C. Comparison with the r-space method 

In order to check the validity of solving the Dirac equa- 
tion in the WS basis and to provide numerical experiences 
for future applications, we compare results of 16 and 
208 Pb from solving the Dirac equation in the WS basis 
and those from solving the same equation in the coordi- 
nate space. The latter is the most accurate method of 
solving the Dirac equation for realistic nuclei up to now 



thus is used as a standard here. The scalar and vector 
potentials in the Dirac equation are provided by very ac- 
curate SRHR calculations with the parameter set NL3 
for the Lagrangian, the mesh size Ar = 0.05 fm, the box 
size i? max = 30 fm for 16 and i? max = 35 fm for 208 Pb. 
Then with thus obtained S(r) and V(r), the Dirac equa- 
tion is solved in the coordinate space and in the WS basis 
with also the parameter set NL3. 

To compare results only from solving the Dirac equa- 
tion avoids errors from other numerical procedures, e.g., 
the error from the iteration and that from solving the 
Laplace equations. For the same reason, what we com- 
pare between these two methods is not the binding energy 
which contains the contribution from mesons but the av- 
erage single particle energy E sp /A. E sp = J2i e i where Ci 
is the single particle energy and the summation over all 
occupied states for both neutrons and protons. We also 
compare the rms radius (r 2 ) 1 / 2 and (r 4 ) 1 / 4 . The radius 
(r 4 ) 1 / 4 reflects the nucleon densities in the tail more than 
the rms radius. 

Table [I] presents the dependence of results of the Dirac 
equation in the SWS basis on the mesh size Ar. With 
Ar decreasing, results in the SWS basis approach the 
standard results, i.e., those in the r space. Ar = 0.1 fm 
gives results accurate enough. The dependence on the 
box size i? max and on the basis size determined by E cut 
are investigated and shown in Figs. ^ and [21 where the 
deviations of the average single particle energy E sp /A, 
the rms radius (r 2 ) 1 / 2 and (r 4 ) 1 / 4 from the standards 
are plotted versus E cut for different i? max - If -Rmax is 
not large enough, it is difficult to approach the standard 
results. For example, when i? max = StqA 1 ^ 3 — 9.4 fm 
and E cut = 300 MeV for 16 O (correspondingly, 7V max ~ 
20), the results seem converge, but the discrepancy of 
the average single particle energy from the standard one 
remains 0.1 keV. So one must use a large enough box with 
the size i? max around 4r A 1,/3 for light nuclei and StqA 1 ^ 3 
for heavy ones. It is interesting that the convergence of 
the results does not depend on N max but only on E cut . 
For 16 ( 208 Pb), the results converge to the standard 
ones at ~ 300 (400) MeV. From Figs. □ and El we find 
the radius (r 4 ) 1 ' 4 also converges very well which implies 
that nucleon densities could be calculated accurately even 
for large r. 

We have made similar investigations for results in the 
DWS basis and similar conclusions are made. For in- 
stances, the deviations of the average single particle en- 
ergy E sp /A, the rms radius (r 2 ) 1//2 and (r 4 ) 1//4 from the 
standards are plotted versus E cut for different i? max in 
Figs. and H 

In the expansion of the nucleon wave function, 
Eq. (|21|l . one has to take into account not only the lev- 
els in the Fermi sea but also those in the Dirac sea be- 
cause they form a complete basis together. Now the ques- 
tion arises how many levels in the Dirac sea one has to 
take into account. In the calculations in Figs. [3] and 01 
we have used n max = « max with rt max determined by 
£" cu t. In Table ITTT1 the dependence of the average sin- 
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TABLE III: Dependence of the single particle energy, rms ra- 
dius, and (r 4 ) 1//4 on the maximum principal quantum number, 
-^max = 2n max + I, for the SRHDWS theory. The meson and 
Coulomb fields are obtained from SRHR calculations with the 
parameter set NL3, Ar = 0.05 fm, _R m ax = 30 fm for le O and 
35 fm for 208 Pb. In SRHDWS calculations, the parameter set 
NL3 is used. For le O, R max = 4r A 1/3 = 12.5 fm, E cut = 300 
MeV for positive states. For 208 Pb, R max = 3r A 1/3 = 22.1 
fm, -Ecut = 200 MeV for positive states. 



N max 


-E sp /A 


(r 2 ) 1 / 2 


( r y/4 


-E sp /A 


(r 2 ) 1 / 2 


(r y/4 






v ) 






2U8 Ph 




no 


23.1129 


2.5912 


2.8859 


23.3331 


5.6314 


5.9889 





23.1077 


2.5916 


2.8861 


23.3329 


5.6314 


5.9889 


2 


23.0762 


2.5939 


2.8889 


23.3316 


5.6315 


5.9890 


4 


23.0617 


2.5942 


2.8893 


23.3304 


5.6317 


5.9892 


6 


23.0439 


2.5946 


2.8898 


23.3299 


5.6318 


5.9893 


8 


23.0385 


2.5946 


2.8899 


23.3294 


5.6319 


5.9893 


10 


23.0376 


2.5946 


2.8899 


23.3292 


5.6319 


5.9894 


12 


23.0375 


2.5946 


2.8899 


23.3291 


5.6319 


5.9899 


14 


23.0375 


2.5946 


2.8899 


23.3290 


5.6319 


5.9894 


16 


23.0375 


2.5946 


2.8899 


23.3290 


5.6319 


5.9894 


18 


23.0375 


2.5946 


2.8899 


23.3289 


5.6319 


5.9893 


20 


23.0375 


2.5946 


2.8899 


23.3288 


5.6319 


5.9893 


22 


23.0375 


2.5946 


2.8899 


23.3287 


5.6319 


5.9893 


30 


23.0375 


2.5946 


2.8899 


23.3287 


5.6319 


5.9893 



gle particle energy, the rms radius (r 2 ) 1 / 2 and (r 4 ) 1 / 4 on 
Af max = 2n^ ax + I — a cutoff on the principal quantum 
number of levels in the Dirac sea — are given for 16 O 
and 208 Pb. From Table ITTTl we find a merit of solving the 
Dirac equation in the DWS basis: the number of nega- 
tive states included in the basis could be much smaller 
that that of the positive states. Let take 16 O as an ex- 
ample, i?max = 4r A 1 / 3 and E cut = 300 MeV for positive 
states correspond to iV+ ax ~ 28. But for negative states, 
Af max = 10 gives very accurate results, e.g., the discrep- 
ancy of E sp /A from the standard is smaller than 0.1 keV. 
This will significantly simplify the deformed problem by 
decreasing the matrix dimension compared to solve the 
Dirac equation in the SWS basis. 

The above investigations are somehow academic. In 
practical applications, it is not necessary to go to the 
accuracy around keV in the single particle energy or 10 
fm in the radius. So in the following calculations, we 
will use i? max = 20 fm, Ar = 0.1 fm and E cut = 60-80 
MeV for heavy and light nuclei which give reasonable 
accuracies both for the binding energy and the radius. 
This set of cutoff's corresponds approximately iV max = 
2n max + 1 ~ 25 where I is the orbital angular momentum 
of relevant state. 



IV. RESULTS AND DISCUSSIONS 

In this section we present results of SRHWS. Since 
our main aim is to show the virtues of SRHWS com- 
pared to SRHHO and SRHR, we do not include pairing 



correlations and restrict our study to doubly magic or 
magic nuclei only. If not specified, the parameter set 
NLSH is used for the Lagrangian, i? max = 20 fm and 
Ar = 0.1 fm throughout this section. Other parame- 
ter sets for the Lagrangian do not change the conclu- 
sion here. In SRHDWS, the number of normal levels in 
the Fermi sea and that of negative ones in the Dirac sea 
are the same for convenience, i.e., rt+ ax = n~ ax . For 
SRHHO, hu = 41 J 4~ 1 / 3 has been used and cutoff's of 
the expansion for fermions and bosons are the same, i.e., 
N F = N B = AW 



A. Bulk properties of stable nuclei from different 
SRH theories 

In Table ITVl the binding energy per nucleon (E/A) and 
neutron, proton and charge radii (r n , r p and r c ) of some 
typical spherical nuclei are presented which are calcu- 
lated from the present available codes, including SRHR 
SRHSWS, SRHDWS, SRHHO. Available data [H Hf 
are also included for comparison. We use approximately 
the same AT max in the SRHHO as that in the SRHWS 
which is determined by E cut . 

Generally speaking, for each studied nucleus, the four 
approaches give almost the same results with an accuracy 
within 0.1% with few exceptions where the differences are 
still less than 0.3%. They are in excellent agreement with 
available data. 

With the same parameters of spatial lattice i? max and 
Ar, SRHWS should reproduce results of SRHR when 
-Ecut (or N max ) is large enough. This is justified in 
Table IIVI One find exactly coincident results between 
SRHSWS and SRHR for most of the studied nuclei. The 
remaining differences and those between SRHDWS and 
SRHR could be diminished by increasing E cut . 

In Figures and {7\ the neutron density distributions 
are compared between SRHR, SRHSWS, SRHDWS, and 
SRHHO, in which 16 O, 48 Ca and 208 Pb are chosen as 
examples. The calculation details are the same as Ta- 
ble IIVI For these stable nuclei, all these SRH methods 
are valid and all calculations arc in excellent agreement 
with each other from the central to outer region of each 
nucleus. Tiny differences in the central region do not af- 
fect the physical observables like the binding energy or 
nuclear radius as is seen in Table ITVl Furthermore, these 
differences could also be decreased by increasing E cu % or 

^max • 

From the above discussions, it is clear that SRHWS 
is equivalent to SRHR and SRHHO for stable nuclei. 
Thus we conclude that Woods-Saxon basis provide an- 
other possibility to solve (non-)relativistic mean field the- 
ory. 
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TABLE IV: The binding energy per nucleon and neutron, 
proton and charge radii of some typical spherical nuclei. The 
parameter set NLSH is used for the Lagrangian. 7? max = 20 
fm and Ar = 0.1 fm for SRHR and SRHWS. hio = 41A~ 1/3 
for SRHHO. Numbers in brackets in the second column give 
E cut for SRHWS and iV max for SRHHO. Data for E/A and 
r c are taken from Ref. |2o| and Ref. [gfj, respectively. Energy 
is in MeV and radius in fm. 
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B. Neutron density distributions for Ca in 
different SRH theories 




r [fm] 



FIG. 5: Neutron density distributions for O from different 
SRH approaches. The parameter set NLSH is used for the 
Lagrangian. i? max = 20 fm and Ar = 0.1 fm for SRHR and 
SRHWS. Scut = 80 MeV for SRHWS. Correspondingly, iV max 
= 25 for SRHHU. In SRHDWS, the number of levels in the 
Dirac sea included in each block is the same as that of nor- 
mal levels which is determined by -E C ut- The inset presents 
logarithmic densities. 




/-[fm] 



FIG. 6: Neutron density distributions for Ca from differ- 
ent SRH approaches. The parameter set NLSH is used for 
the Lagrangian. i? max = 20 fm and Ar = 0.1 fm for SRHR 
and SRHWS. E cut = 80 MeV for SRHWS. Correspondingly, 
iV max = 25 for SRHHU. In SRHDWS, the number of levels 
in the Dirac sea included in each block is the same as that of 
normal levels which is determined by _E cut . The inset presents 
logarithmic densities. 



As we already discussed in the introduction, one of the 
merits of SRHR against SRHHO is its proper description 
of exotic nuclei. In this subsection, we will demonstrate 
the equivalence between SRHWS and SRHR when rea- 
sonably large E cut is applied in SRHWS. 

In order to see the results for the unstable nuclei near 



the neutron drip line, the neutron density distribution for 
72 Ca is studied here. The nucleus 72 Ca is predicted to 
be the last bound calcium isotope [23, |28|, |29, |3(j . Since 
it is not a doubly magic nucleus, there might be some 
uncertainty in present results due to the lack of inclusion 
of pairing correlations. However, as we stressed in the 
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FIG. 7: Neutron density distributions for Pb from differ- 
ent SRH approaches. The parameter set NLSH is used for 
the Lagrangian. i? max = 20 fm and Ar = 0.1 fm for SRHR 
and SRHWS. E cut = 60 MeV for SRHWS. Correspondingly, 
iV m ax = 25 for SRHHO. In SRHDWS, the number of levels 
in the Dirac sea included in each block is the same as that of 
normal levels which is determined by E cut . The inset presents 
logarithmic densities. 



TABLE V: Convergence study for 72 Ca. The parameter set 
NLSH is used for the Lagrangian. Ar = 0.1 fm for SRHR and 
SRHWS. Scut = 75 MeV for SRHWS. Energy is in MeV and 
radius in fm. 





E/A 


r n 




An 


^max 






SRHR 




20 


6.482 


4.656 


3.639 


-0.191 


25 


6.483 


4.723 


3.639 


-0.221 


30 


6.484 


4.773 


3.639 


-0.228 


35 


6.484 


4.807 


3.639 


-0.229 


^max 




SRHSWS 




20 


6.481 


4.663 


3.639 


-0.206 


25 


6.482 


4.726 


3.639 


-0.231 


30 


6.483 


4.774 


3.639 


-0.237 


35 


6.483 


4.803 


3.639 


-0.238 


Rmax 




SRHDWS 




20 


6.474 


4.662 


3.641 


-0.163 


25 


6.475 


4.733 


3.641 


-0.197 


30 


6.475 


4.789 


3.640 


-0.205 


35 


6.475 


4.828 


3.640 


-0.206 


iV max 






SRHHO 




25 


6.489 


4.577 


3.639 


-0.054 


31 


6.492 


4.605 


3.639 


-0.128 


37 


6.494 


4.628 


3.639 


-0.166 


43 


6.494 


4.649 


3.639 


-0.189 



beginning of this section, the main aim here is to show 
the virtue of SRHWS compared to SRHHO, it is very 
unlikely that the pairing would change our conclusion 
qualitatively. 

For stable nuclei, it has been shown that i? max ~ 20 
fm is large enough. For drip line nuclei, the dependence 



FIG. 8: Comparison of density distributions for Ca from 
SRHR, SRHSWS and SRHHO. The parameter set NLSH 
is used for the Lagrangian. Ar = 0.1 fm for SRHR and 
SRHSWS. E cut = 75 MeV and i? max = 20 (thick dot-dashed 
curve), 25 (thick long-dashed curve), 30 (thick dashed curve) 
and 35 fm (thick dotted curve) for SRHSWS. These sets of 
cutoff's correspond to cutoff's in principal quantum number 
A^max = 25 (thin dot-dashed curve), 31 (thin long-dashed 
curve), 37 (thin dashed curve) and 43 (thin dotted curve) 
which are used in SRHHO calculations. The density distribu- 
tion from SRHR are almost identical with that from SRHSWS 
with the same box size. For brevity, only p n (r) from SRHR 
with fl ma x = 35 fm (thick solid line) is displayed here which 
covers the curve corresponding to p n (r) from SRHSWS with 
-Rmax = 35 fm (thick dotted curve). 



of the results on i? max for 72 Ca is presented in Table [V] 
For both SRHR and SRHWS, Ar = 0.1 fm and i? max 
= 20, 25, 30 and 35 fm have been adopted respectively. 
The energy cutoff E cut = 75 MeV is used to SRHWS. In 
SRHR and SRHWS calculations, the neutron rms radius 
r n and Fermi energy A n of 72 Ca converge around i? max 
= 35 fm while a independence of the binding energy per 
nucleon E/A and proton rms radius r p on the box size 
can be seen. These sets of parameters in SRHWS, E cut 
— 75 MeV and i? max = 20, 25, 30 and 35 fm, correspond 
to cutoff's on principal quantum number N max — 25, 31, 
37 and 43 which are used in SRHHO calculations in order 
to make fair comparisons between SRHWS and SRHHO. 
Similar as those from SRHR and SRHWS, E/A and r p 
depend little on N max in SRHHO. However, the neutron 
rms radius r n increases with N max steadily, which shows 
a much slower convergence. As it is based on a complete 
basis, SRHHO can reach convergence of r n if N max is 
large enough. From Table El one finds that for the same 
N max (or equivalent i? max ), a difference of Ar n 0.2 
fm between SRHHO and SRHWS (SRHR) can be seen. 
From the slow convergence of r n with N max in SRHHO 
(AiV max = 6 gives Ar n w 0.02 fm), we can estimate the 
lower limit of the iV max as iV max « 90 in order to give r n 
= 4.8 fm in SRHR or SRHWS. 
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FIG. 9: Comparison of density distributions for 72 Ca from 
SRHR, SRHDWS and SRHHO. The parameter set NLSH 
is used for the Lagrangian. Ar = 0.1 fm for SRHR and 
SRHDWS. E cut = 75 MeV and i? max = 20 (thick dot-dashed 
curve), 25 (thick long-dashed curve), 30 (thick dashed curve) 
and 35 fm (thick dotted curve) for SRHDWS. In SRHHO 
calculations, iV max = 25 (thin dot-dashed curve), 31 (thin 
long-dashed curve), 37 (thin dashed curve) and 43 (thin dot- 
ted curve) are used. The density distribution from SRHR are 
very close to that from SRHSWS with the same box size. For 
brevity, only p n (r) from SRHR with i? max = 35 fm (thick 
solid line) is displayed here. 



TABLE VI: Effects of negative levels on bulk properties in 
SRHDWS for 16 O. The parameter set for the Lagrangian is 



We compare the neutron density distribution of 72 Ca 
from different SRH approaches in Figs. |<3 and H3 With 
the same box size, the density distribution from SRHR 
are almost identical with those from SRHWS, which in- 
dicates the equivalence between SRHWS and SRHR. For 
brevity, only p n (r) from SRHR with i? max = 35 fm is 
displayed in Figs. |S] or [5] which covers the curve corre- 
sponding to p n {r) from SRHSWS with i? m ax = 35 fm 
in Fig. EI On the other hand, p n {r) from SRHHO even 
with iV max = 43 fails to reproduce the result of SRHR 
due to the localization property of HO potential. In addi- 
tion, with the same -/V max , the spatial extension of p n {r) 
from SRHWS are always larger than that from SRHHO. 
The variational tendency of the curve p n (r ) ~ r also ex- 
plains different convergence behaviors of r n in SRHWS 
and SRHHO as given in Table Ivl With increasing -R max , 
the p n (r) in SRHWS has the correct asymptotic behavior 
while that in SRHHO decay too quickly. 



This result is very encouraging and tell us that even 
the long tail (or halo) behavior in neutron density dis- 
tribution for nuclei near the drip line can be described 
by SRHSWS as well as SRHR, if pairing correlations is 
incorporated properly. 



NLSH, R„ 



20 fm, Ar = 0.1 fm, and N, 



+ 



25. For the 



initial Woods-Saxon like potentials, parameters in Ref. |23f 
are used except for Vb which is specified in the table. The left 
value in each entry gives the result without negative levels 
included and the right one that with iV~ ax = 25. Energy is 
in MeV and radius in fm. 





Vb = 54 MeV 


Vb = 72 MeV 


Vb = 90 MeV 


E/A 


8.547 | 8.013 


8.117 | 8.015 


8.427 | 8.012 


^rms 


2.385 2.568 


2.531 | 2.567 


2.610 2.567 



Contribution from negative levels in the 
SRHDWS theory 



In the expansion of the nucleon wave function, 
Eq. (|21|1 . one has to take into account not only the levels 
in the Fermi sea but also those in the Dirac sea because 
they form a complete basis together. We study the ef- 
fects of negative levels and give results of an example, 
16 O, in Table I VII Firstly, without negative levels in- 
cluded, the nucleus is over bound and the nuclear size 
is smaller as seen from Table IVII Secondly, in contrary 
to the case with negative levels included, the calculated 
nuclear properties depend on the initial potentials very 
much if no negative levels included. 

It should be noted that the contribution from negative 
levels depends on the initial Woods-Saxon potentials for 
generating the DWS basis. So do the cutoff N max or E cnt 
for convergence. If the initial Woods-Saxon potential is 
exactly identical to the converged potentials, the matrix 
in Eq. (12211 is diagonal, negative states do not contribute 
because of the no sea approximation. Positive states can 
also be chosen as less as possible, e.g., lsi/ 2 , IP3/2 and 
lPi/2 are enough for 16 O. From the third column corre- 
sponding to Vb = 72 MeV in Table |Vl| one finds that the 
initial nuclear potential for the Dirac equation proposed 
in Ref. [23| is a good choice for SRHDWS as the negative 
states only contribute ~ 1.25 % to both E/A and r rms . 
If we change the initial potentials, e.g., by changing Vq 
by 25%, much larger contributions from negative states 
are found in Table l*VTl 



In order to know the contribution of negative levels in 
the Dirac sea to the wave function, the value of J 2n \ c n 1 2 
in the expansion, Eq. I|21|) has been given in Table IVlTl for 
occupied states of 16 O. We note that the nucleon wave 
function is normalized to one. It can be seen that a small 
component from negative states in the wave functions, 
about 1CP 4 5 , contributes to the physical observables 
such as E/A and r rms by the magnitude of 1%~10% as 
seen from Table IVII Again we notice that the initial 
Woods-Saxon potentials differ more from the converged 
ones, the larger is the contribution from negative levels. 
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TABLE VII: Contribution of negative levels in the Dirac sea 
to the single nucleon wave functions of 16 O in SRHDWS. For 
each single particle levels, \ c n 1 2 m the expansion Eq. 12111 
is presented. The parameter set for the Lagrangian is NLSH, 
-Rmax = 20 fm and Ar = 0.1 fm and iV+ ax = iV~ ax = 25. 
For the initial Woods-Saxon like potentials, the parameters 
in Ref. |23|| are used except for Vo which is specified in the 
table. The binding energy per nucleon and nuclear radius for 
each calculation is presented in Table IVT1 



Level 


Vo = 54 MeV 


V = 72 MeV 


Vo = 90 MeV 


NLs 1/2 


9.15xl0~ 5 


3.48xl0~ 5 


3.44x10" 


-4 


Nlp 3/2 


5.32xl0 -5 


6.05xl0" 5 


3.69x10" 


-4 


Nlpi/ 2 


4.01xl0~ 4 


1.06xl0~ 5 


8.80x10" 


-4 


PlSl/2 


1.92xl0 -4 


4.92xl0 -5 


4.89x10" 


-5 


P1P3/2 


1.17xl0~ 4 


3.27xl0~ 5 


8.07x10" 


-5 


Plpi/2 


8.09xlO" 4 


2.08xl0~ 4 


2.38x10" 


-5 




V. 


SUMMARY 







We have solved the spherical relativistic Hartree the- 
ory in the Woods-Saxon basis (SRHWS). The Woods- 
Saxon basis is obtained by solving either the Schrodingcr 
equation (SRHSWS) or the Dirac equation (SRHDWS). 
Formalism and numerical details for both cases are pre- 
sented. The WS basis in the SRHDWS theory could be 
very smaller than that in the SRHSWS theory. This will 
largely facilitate solving the deformed problem. 

The results from SRHWS are compared with those 
from solving the spherical relativistic Hartree theory in 
the harmonic oscillator basis, SRHHO, and those in the 
coordinate space, SRHR. For stable nuclei, all approaches 
give identical results for properties such as total binding 
energies and the neutron, proton and charge rms radii as 
well as neutron density distributions. 

For neutron drip line nuclei, e.g. 72 Ca, which has a 
very small neutron Fermi energy A n ~ 0.2 MeV, both 



SRHR and SRHWS easily approach convergence by in- 
creasing box size while SRHHO does not. Furthermore, 
SRHWS can satisfactorily reproduce the neutron density 
distribution from SRHR, but SRHHO fails with similar 
cutoff's. 

In SRHDWS calculations, the negative levels in the 
Dirac Sea must be included in the basis in terms of 
which nucleon wave functions are expanded. We stud- 
ied in detail the effects and contributions of negative 
states. Without negative levels included, the calculated 
nuclear properties depend on the initial potentials very 
much. A small component from negative states in the 
wave functions, about 10 -4 5 , contributes to the phys- 
ical observables such as E/A and r rms by the magnitude 
of 1%~10%. When the initial potentials differ more from 
the converged ones, the contribution from negative levels 
becomes more important. 

We conclude that the Woods-Saxon basis provides a 
reconciler between the harmonic oscillator basis and the 
coordinate space which may be used to describe exotic 
nuclei both properly and efficiently. 

The extension of spherical relativistic Hartree theory 
in the Woods-Saxon basis to deformed cases with pairing 
included is in progress. 
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